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Abstract 

In  the  context  of  methodologies  intended  to  confer  robustness  to^  geometric^  algo¬ 
rithms,  we  elaborate  on  the  exact  computation  paradigm  and  formalize  the  notion  of 
degree  of  a  geometric  algorithm,  as  a  worst-case  quantification  of  the  precision  (num¬ 
ber  of  bits)  to  which  arithmetic  calculation  have  to  be  executed  in  order  to  parantee 
topological  correctness.  We  also  propose  a  formalism  for  the  expeditious  evaluation  of 
algorithmic  degree.  As  an  application  of  this  paradigm  and  an  illustration  of  our  general 
approach,  we  consider  the  important  classical  problem  of  proidmity  queries  in  2  and  3 
dimensions,  and  develop  a  new  technique  for  the  efficient  and  robust  execution  of  sue 
queries  based  on  an  implicit  representation  of  Voronoi  diagrams.  Our  new  technique 
gives  both  low  degree  and  fast  query  time,  and  for  2D  queries  is  optimal  with  respect  to 
both  cost  measures  of  the  paradigm,  asymptotic  number  of  operations  and  arithmetic 
degree. 
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1  Introduction 

The  increasing  demand  for  efficient  and  reliable  geometric  software  libraries  in  key  appli¬ 
cations  such  as  computer  graphics,  geographic  information  systems,  and  computer-aided 
manufacturing  is  stimulating  a  major  renovation  in  the  field  of  computational  geometry. 
The  inadequacy  of  the  traditional  simplified  framework  has  become  apparent,  and  it  is 
being  realized  that,  in  order  to  achieve  an  effective  technology  transfer,  new  frameworks 
and  models  are  needed  to  design  and  analyze  geometric  algorithms  that  are  efficient  in  a 
practical  realm. 

The  real-RAM  model  with  its  implicit  infinite-precision  requirement,  has  proved  un¬ 
realistic  and  needs  to  be  replaced  with  a  realistic  finite-precision  model  where  geometric 
computations  can  be  carried  out  either  exactly  or  with  a  guaranteed  error  bound.  This  has 
motivated  a  great  deal  of  research  on  the  subject  of  robust  computational  geometry  (see, 
e.g.,  [4,  11,  10,  17,  23,  24,  26,  29,  28,  32,  39,  45,  48,  19]).  Also,  efficiency  must  be  evalu¬ 
ated  in  a  finer  framework  than  the  conventional  “big-Oh”  analysis.  In  particular,  constant 
factors  dependent  on  the  precision  requirement  of  the  numerical  computations  should  be 
taken  into  account.  For  an  early  survey  of  the  different  approaches  to  robust  computational 
geometry  the  reader  is  referred  to  [31]. 

To  a  first,  rough,  approximation,  robustness  approaches  are  of  two  main  types:  per¬ 
turbing  and  nonperturbing.  Perturbing  approaches  transform  the  given  problem  into  one 
that  is  intended  not  to  suffer  from  well-identified  shortcomings;  nonperturbing  approaches 
are  based  on  the  notion  of  “exact”  (rather  than  “approximate”)  computations,  with  the 
assumption  that  (bounded-length)  input  data  are  error  free.  In  this  category  falls  the 
exact  geometric  computation  paradigm  independently  advocated  by  Yap  [49]  and  by  the 
Saarbriicken  school  [9],  and  so  does  our  approach.  Within  this  paradigm,  we  introduce 
the  concept  of  degree  of  a  geometric  algorithm,  which  characterizes,  up  to  a  small  additive 
constant,  the  arithmetic  precision  (i.e.,  number  of  bits)  needed  by  a  large  class  of  geo¬ 
metric  algorithms.  Namely,  if  the  coordinates  of  the  input  points  of  a  degree-d  algorithm 
within  this  class  are  6-bit  integers,  then  the  algorithm  may  be  required  on  some  instances 
to  perform  arithmetic  computations  with  bit  precision  d6-|-0(l). 

Theoretical  analysis  and  experimental  results  show  that  multiprecision  numerical  com¬ 
putations  take  up  most  of  the  CPU  time  of  exact  geometric  algorithms  (see,  e.g.,  [34,  40]). 
Thus,  we  believe  that,  in  defining  the  efficiency  of  a  geometric  algorithm,  the  degree  should 
be  considered  as  important  as  the  asymptotic  time  complexity  and  should  correspondingly 
play  a  major  role  in  the  design  stage. 

Research  along  these  lines  involves  re-examining  the  entire  rich  body  of  computational 
geometry  as  we  know  it  today.  In  this  paper,  we  consider  as  a  test  case  a  problem  area, 
geometric  proximity,  which  plays  a  major  role  in  applications  and  has  recently  attracted 
considerable  attention  because,  due  to  its  demands  of  high  precision  for  exact  computation, 
it  is  particularly  appropriate  in  assessing  effectiveness  of  robust  approaches  (see,  e.g.,  [8, 
10,  18,  25,  26,  23,  46]). 

While  Voronoi  diagrams  are  interesting  in  their  own  right,  the  main  reason  for  con- 
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strutting  and  storing  them  is  to  efficiently  answer  fundamental  proximity  queries  such  as 
nearest-neighbor  and  circular-range  queries.  In  this  paper,  we  use  the  notion  of  degree  to 
illustrate  the  drawbacks  of  current  approaches,  and  adopt  it  as  a  design  criterion  in  deve  - 
oping  a  new  technique  for  the  efficient  and  robust  execution  of  proximity  queries,  based  on 

an  implicit  representation  of  Voronoi  diagrams. 

Throughout  this  paper,  we  assume  that  the  coordinates  of  all  input  data  (also  called 
primitive  points)  are  f»-bit  integers.  Hence,  parameter  b  denotes  the  input  precision.  In  order 
to  perform  exact  computations,  the  coordinates  of  the  points  computed  by  the  algorithm 
(referred  to  here  as  derived  points,  e.g.,  the  vertices  of  a  Voronoi  diagram  of  points  and 
segments)  must  be  stored  with  a  representation  scheme  that  supports  rational  or  algebraic 
numbers  as  data  types  (through  multiprecision  integers). 

Consider,  for  example,  the  nearest-neighbor  query  problem:  given  a  set  5  of  n  sites  in  the 
plane,  find  the  site  closest  to  a  query  point  q.  The  well-known  optimal  deterministic  method 
for  this  problem  consists  of  performing  point-location  in  the  Voronoi  diagram  of  5,  denoted 
V{S).  This  takes  optimal  time  O(logn).  However,  in  the  straightforward  implementation 
of  this  method,  the  coordinates  (x,  y)  of  a  vertex  of  V (5)  are  rational  numbers  given  by  the 
ratio  of  two  3x3  determinants  whose  entries  are  integers  of  well-defined  maximum  modulus. 
Hence,  homogeneous  coordinates  iX,Y,W)  are  used  for  the  exact  representation  of  y(5), 
where  i  =  X/W,  y  =  Y/W,  X  and  Y  are  3&-b  bit  integers,  and  W  is  a  26-bit  integer.  The 
fundamental  operation  used  by  any  point  location  algorithm  is  the  point-line  discrimination, 
which  consists  of  determining  whether  the  query  point  q  is  to  the  left  or  to  the  right  of  an 
edge  between  vertices  ui  and  V2-  For  the  case  of  the  Voronoi  diagram  V’(S),  this  is  equivalent 
to  evaluating  the  sign  of  a  3  X  3  determinant  whose  rows  are  the  homogeneous  coordinates  of 
q.  vu  and  U2,  a  computation  that  needs  about  66  bits  of  precision.  This  should  be  compared 
with  the  0(7i)-time  brute-force  method  that  computes  the  (squares  of  the)  distances  from 
q  to  all  the  sites  of  S,  and  executes  arithmetic  computations  with  only  26  bits  of  precision 
(which  is  optimal). 

The  technique  presented  in  this  paper  uses  a  new  point-location  data  structure  for 
Voronoi  diagrams  such  that  the  test  operations  executed  in  the  point  location  procedure 
are  executed  with  optimal  26-bit  precision.  Hence,  our  approach  reconciles  efficiency  with 

robustness. 

Our  approach  also  supports  an  object-oriented  programming  style  where  access  to  the 
geometry  of  Voronoi  diagrams  in  point-location  queries  is  encapsulated  in  a  small  set  of 
geometric  test  primitives. 

The  main  results  of  this  work  are  summarized  in  Table  1.  Considering,  for  the  time 
being,  the  degree  as  a  measure  of  complexity,  we  show  that  previous  methods  exhibit  a 
sharp’ tradeoff  between  degree  and  query  time.  Namely,  low  degree  is  achieved  by  the  slow 
brute-force  search  method,  while  fast  algorithms  based  on  point-location  in  a  preprocessed 
Voronoi  diagram  have  high  degree.  Our  new  technique  gives  instead  both  low  degree  and 

’While  this  research  was  nearing  completion,  Bumikel’s  thesis  [8]  was  distributed  on  the  Web,  and  we 
became  aware  that  Bumikel  had  independently  defined  the  concept  of  precision,  which  is  equivalent  of  our 
notion  of  degree. 
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fast  query  time,  and  is  optimal  with  respect  to  both  cost  measures  for  2D  queries  in  sets  of 
points. 


Query 

Method 

Degree 

Time 

nearest  neighbor 

brute-force  distance  comparison 

2  ♦ 

0(n) 

point  location  in  explicit  Voronoi  diagram 

6 

0(lognj  * 

point  location  in  implicit  Voronoi  diagram 

2  * 

M  IIMlilMII 

/c-nearest  neighbors 
and 

circular  range  search 

brute-force  distance  comparison 

2* 

0(n) 

point  location  in  explicit  order-fc  Voronoi  diagram 

6 

point  location  in  implicit  order-k  Voronoi  diagram 

2  * 

0(logn  -f  k)  *  II 

nearest  neighbor  among 
points  and  segments 

brute-force  distance  comparison 

6 

point  location  in  explicit  Voronoi  diagram 

64 

■BfiESSlill 

point  location  in  implicit  Voronoi  diagram 

6 

3D  nearest  neighbor 

brute-force  distance  comparison 

2  ♦ 

0(n) 

point  location  in  explicit  3D  Voronoi  diagram 

8 

0{\og^  n) 

point  location  in  implicit  3D  Voronoi  diagram 

3 

0(log^ 

Table  1:  Comparison  of  the  degree  and  time  of  algorithms  for  some  fundamental  proximity 
query  problems.  A  *  denotes  optimality.  The  new  technique  introduced  in  this  paper 
{point  location  in  an  implicit  Voronoi  diagram)  always  outperforms  previous  methods  and 
is  optimal  for  2D  queries. 

The  rest  of  this  paper  is  organized  as  follows.  Our  approach  is  described  in  Section  2, 
where  the  concept  of  degree  of  a  geometric  algorithm  is  formalized.  In  Section  3,  we  consider 
the  following  fundamental  proximity  queries  for  a  set  of  point  sites  in  the  plane;  nearest 
neighbor  search,  fc-nearest  neighbors  search,  and  circular  range  search.  We  show  that  the 
existing  methods  for  answ'ering  such  queries  efficiently  have  degree  6,  and  we  present  our 
new  technique,  based  on  implicit  Voronoi  diagrams,  that  achieves  optimal  degree  2.  In 
Sections  4-5,  we  extend  our  approach  to  nearest  neighbor  search  queries  in  a  set  of  3D 
point  sites  and  in  a  set  of  point  and  segment  sites  in  the  plane,  respectively.  Practical 
improvements  are  presented  in  Section  6.  Finally,  further  research  directions  are  discussed 
in  Section  7. 

2  Degree  of  Geometric  Problems 

The  numerical  computations  of  a  geometric  algorithm  are  basically  of  two  types;  tests 
and  constructions.  The  two  types  of  computations  have  clearly  distinct  roles.  Tests  are 
associated  with  branching  decisions  in  the  algorithm  that  determine  the  flow  of  control, 
whereas  constructions  are  needed  to  produce  the  output  data  of  the  algorithm. 

Input  data  (whether  the  result  of  empirical  observations  or  not)  are  reasonably  assumed 
to  be  expressed  with  h  bits,  for  some  small  integer  h.  Approximations  in  the  execution  of 
constructions  give  rise  to  approximate  results,  which  may  be  not  only  acceptable  but  also 
mandated;  indeed,  as  long  as  their  maximum  absolute  error  does  not  exceed  the  resolution 
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required  by  the  application  (such  as  spacing  of  raster  lines  in  graphics),  it  would  be  wasteful 

to  produce  substantially  more  accurate  results. 

On  the  other  hand,  approximations  in  the  execution  of  tests  may  produce  an  incorrect 
branching  of  the  algorithm.  Such  event  may  have  catastrophic  consequences,  giving  rise  to 
structurally  incorrect  results.  Therefore,  tests  are  much  more  critical,  and  their  execution 

must  be  carried  out  with  complete  accuracy. 

We  shall  therefore  characterize  geometric  algorithms  on  the  basis  of  the  complexity  of 
their  test  computations.  In  a  large  class  of  geometric  algorithms,  tests  are  based  on  the 
evaluation  of  the  sign  of  a  multivariate  polynomial.  Notice  that,  in  general,  when  all  input 
data  are  dimensionally  homogeneous  (such  as  coordinates  of  points),  the  polynomial  is  also 
homogeneous  because  each  of  its  monomials  has  the  physical  dimension  of  a  voluine  in  d- 
dimensional  space  (for  some  integer  d  depending  upon  the  nature  of  the  test).  Expeninental 
results  (see,  e.g.,  [34,  40])  show  that  a  substantial  fraction  of  the  CPU  time  is  spent  in  the 
evaluation  of  such  polynomials. 

Now  we  define  the  concept  of  arithmetic  degree.  We  consider  algorithms  that  evaluate 
multivariate  polynomials  over  their  variables.  A  primitive  variable  is  an  input  variable  of  the 
algorithm  and  has  conventional  arithmetic  degree  1.  The  arithmetic  degree  of  a  variable  u 
is  the  arithmetic  degree  of  the  multivariate  polynomial  E  that  computes  v.  The  arithmetic 
degree  of  E  is  the  maximum  (in  the  homogeneous  case,  the  common)  arithmetic  degree  of 
its  monomials.  The  arithmetic  degree  of  a  monomial  is  the  sum  of  the  arithmetic  degrees  of 
of  its  variables.  As  observed  above,  our  definition  of  arithmetic  degree  coincides  with  that 
of  precision  in  Burnikel’s  work  [8]  and  is  related  to  that  of  depth  proposed  by  Yap  [48,  49], 
as  discussed  below. 

Definition  1  An  algorithm  has  degree  d  if  its  test  computations  involve  the  evaluation  of 
multivariate  polynomials  of  arithmetic  degree  at  most  d. 

We  make  the  assumption  that  every  multivariate  polynomial  used  in  tests  depends 
upon  a  bounded-size  set  of  primitive  variables  and  therefore  has  0(1)  monomials.  This 
assumption  holds  for  a  large  class  of  geometric  computations,  including  those  discussed  in 

this  paper. 

An  immediate  consequence  of  Definition  1  and  of  the  above  assumption  is  the  following 
fact,  which  justifies  our  use  of  the  degree  of  an  algorithm  to  characterize  the  precision 

required  in  test  computations. 

Lemma  1  If  an  algorithm  has  degree  d  and  its  input  variables  are  b-bit  integers,  then  all 
the  test  computations  can  be  carried  out  with  db  -|-  0(1)  bits. 

Definition  2  A  problem  11  has  degree  d  if  any  algorithm  that  solves  11  has  degree  at  least 
d. 


It  is  worth  mentioning  that  related  definitions  have  been  given  by  Yap  and  are  based  on 
the  concept  depth  of  derivation  [48,  49].  Given  a  set  of  numbers,  any  number  x  of  the  set 
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hcis  depth  0.  A  number  has  depth  at  most  d  if  it  can  be  obtained  by  executing  a  rational 
operation  on  numbers  with  depth  d  —  1  or  it  is  the  result  of  a  root  extraction  from  a  degree 
k  polynomial  whose  coefficients  have  depth  at  most  d  —  k.  An  algorithm  has  depth  d  if  it 
performs  only  rational  operations  such  that  all  the  intermediate  computed  numbers  have 
depth  of  derivation  at  most  d  with  respect  to  the  set  of  input  numbers.  Clearly  d  is  the  least 
possible  integer  such  that  all  the  intermediate  computed  values  have  depth  of  derivation  at 
most  d.  A  problem  has  depth  d  if  it  can  be  solved  by  an  algorithm  with  rational  bounded 
depth  d. 

Although  the  definition  of  depth  of  a  problem  can  resemble  the  one  of  degree  of  a 
problem,  the  latter  seems  to  be  more  appropriate  to  the  type  of  study  that  we  want  to 
develop,  where  we  aim  at  minimizing  the  number  of  bits  needed  for  computing  an  exact 
value,  independently  of  its  (possibly  very  high)  depth. 

2.1  Degree  of  an  Algorithm 

Typically  the  support  of  a  geometric  test  is  not  naturally  expressed  by  a  multivariate  polyno¬ 
mial,  but,  rather,  by  a  pair  (£^i,  £'2)  of  expressions  involving  the  four  arithmetic  operations, 
powering,  and  the  extraction  of  square  roots.  The  test  is  typically  the  magnitude  com¬ 
parison  of  £1  and  £2.  Any  such  expression  can  be  viewed  as  a  binary  tree,  whose  leaves 
represent  input  variables  and  whose  internal  nodes  are  of  two  types:  binary  nodes,  which 
are  labeled  with  an  operation  from  the  set  ,  X,  and  unary  nodes,  which  are  labeled 
either  with  a  power  or  with  a  square  root  extraction  (notice  hat  we  restrict  ourselves  to  this 
type  of  radicals).  If  no  radical  appears  in  the  trees  of  £1  and  £2,  then  the  test  is  trivially 
equivalent  to  the  evaluation  of  the  sign  of  a  polynomial,  since  £,*  is  a  rational  function  of 
the  form  ^  {i  =  1,2,  iV,-,  A  are  not  both  trivial  polynomials)  and 


A  >  £2  4=^  (-1)^(^0+<-P2)(aD2  -  N2D,)  >  0 


where 


■’(£)= 


{ 


1  if  E  <  0 
0  if  >  0 


(Note  that  the  above  predicate  implies  the  inductive  assumption  that  the  signs  of  lower 
degree  expressions  Ni,N2,  Di,  and  D2  are  known.) 

Suppose  now  that  at  least  one  of  the  trees  of  Ei  and  E2  contains  radicals.  We  prune 
the  trees  so  that  the  pruned  trees  contain  no  radicals  except  at  their  leaves.  Then  Ni  and 
Di  (i  =  1,2)  can  be  viewed  as  polynomials  whose  variables  are  the  radicals  and  whose 
coefficients  are  (polynomial)  functions  of  non-radicals.  Given  a  polynomial  E  in  &  set  of 
radicals,  for  any  radical  R  in  this  set,  we  can  express  E  as 


E  =  E"R  +  E' 


where  neither  E"  nor  E'  contains  R.  Then 
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E>0  E"R  >  -E'. 

If  E'  and  E2  have  the  same  sign,  then  the  sign  of  E  is  their  common  sign.  Otherwise, 


E  >  a  ^  -  E'^)  >  0 

The  resulting  expression  {E"^R'^  -  E''^)  does  not  contain  R.  Therefore  by  this  device, 
referred  to  as  segregate  and  square,  we  can  eliminate  one  radical.  This  shows  that  by  the 
four  rational  operations  we  can  reduce  the  sign  test  to  the  computation  of  the  sign  of  a 
multivariate  polynomial. 

We  now  present  a  simple  formalism  that  enables  us  to  rapidly  obtain  the  degree  of  the 
sought  multivariate  polynomial  equivalent  to  the  original  algebraic  expression. 

The  terms  involved  in  the  formal  manipulations  are  of  two  types,  generic  and  specific. 
Generic  terms  have  the  form  a*  (for  a  formal  variable  a  and  an  integer  s),  representing 
an  unspecified  multivariate  polynomial  of  degree  $  over  primitive  variables.  Specific  terms 
have  the  form  pj,  for  some  integer  index  j,  representing  a  specified  expression  involving 
the  operators  x,^,vT.  The  terms  are  members  of  a  free  commutative  semiring, 

i.e.,  addition  and  multiplication  are  associative  and  commutative,  addition  distributes  over 
multiplication,  and  specific  terms  can  be  factored  out.  Beside  these  conventional  algebraic 
rules,  we  have  a  set  of  rewriting  rules  of  the  form  A  — ^  B,  meaning  that  the  sign  of  A  is 
uniquely  determined  by  the  sign  of  B.  Note,  the  two  signs  not  always  coincide:  indeed  the 
correspondence  between  the  signs  of  the  sides  depends  upon  the  evaluated  signs  of  other 
expressions  of  lower  degree  (see  the  example  at  the  beginning  of  this  section) . 

The  first  of  these  rules  is  genericization,  i.e.,  a  specific  term  pj,  which  is  known  to  be  a 
polynomial  of  degree  s  over  primitive  variables,  can  be  rewritten  as 

(1)  Pj  — >  a"  (genericization) 

Next,  since  in  this  context  the  only  relevant  feature  of  a  polynomial  is  its  degree,  we 
have 


(2) 

— »  a»+’‘ 

(product  rule) 

(3)  + 

— a* 

(sum  rule) 

(4)  -a* 

—4  a* 

(sign  rule). 

The  role  of  specific  terms  is  that  we  wish  to  keep  track  of  their  structure  (that  is,  their 
definition),  in  order  to  exploit  it  when  computing  least  common  multiples  or  multiplying 
radicals  together.  We  have  therefore  the  following  additional  rewriting  rules 

(5)  fi-  i  ^  ^  Pj  i  Ph  (common  denominator) 

(6)  ^  ^  ^  PjPk  i  PiPh  (common  denominator) 

(7)  Pi±Pj  — >  p]~p]  (segregate  and  square). 
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Note  that  rule  (5)  is  correct  only  in  our  context  of  analysis  of  the  degree  of  geometic 
tests. 

To  illustrate  this  approach  we  discuss  the  simple  example  of  the  point-to-lines  distance 
test,  i.e.,  given  two  lines  rj  and  r2  on  the  plane  and  a  query  point  5,  determine  whether  g 
is  closer  to  rj  than  to  r2. 

Lemma  2  The  point-to-lines  distance  test  can  be  solved  with  degree  6. 

Proof:  Let  the  equation  of  r,  be  a, a:  -t-  biy  ■+■  c,  =  0  (i  =  1, 2)  and  let  g  =  (x,,  yq).  Then  the 
test  is  to  study  the  sign  of  .  By  using  the  proposed  notation, 

and  with  obvious  meaning  for  pi  and  p2,  this  test  becomes  (each  arrow  being  superscripted 
with  the  rules  used) 

_^(6)  Q,2^2-aVi  — — ^(^-3) 

□ 

The  following  lemmas  describe  the  degree  of  other  proximity  primitives  that  will  be 
useful  in  the  rest  of  the  paper.  We  omit  the  proofs  of  such  lemmas,  since  they  are  either 
straightforward  or  have  been  already  proved  in  [8].  However,  it  is  worth  mentioning  that  the 
proofs  in  [8]  can  be  substantially  simplified  by  using  the  notation  proposed  in  this  paper. 

Let  p  be  a  point  and  r  a  line  in  the  plane.  The  point-to-point-line  distance  test  determines 
whether  a  query  point  g  is  closer  to  p  or  to  r. 

Lemma  3  The  point-to-point-line  distance  test  can  be  solved  with  degree  4. 

Let  Pi  and  p2  be  two  distinct  points  of  the  plane  and  let  g  be  a  query  point.  The  points 
distance  test  determines  whether  g  is  closer  to  pi  or  to  p2- 

Lemma  4  The  point-to-points  distance  test  can  be  solved  with  degree  2. 

Notice  that  the  above  lemma  can  be  easily  extended  to  any  space  of  dimension  d. 
Another  fundamental  proximity  primitive  is  the  incircle  test,  that  is  testing  whether  the 
circle  determined  by  three  distinct  sites  (points  and  segments)  of  the  plane  contains  a  given 
query  site.  The  incircle  test  is  a  basic  operation  for  many  algorithms  that  construct  the 
Voronoi  diagram  of  the  sites  (see,  e.g.  [30,  34,  27,  3]).  The  degree  of  the  incircle  test  has 
been  extensively  studied  by  Burnikel  [8]  and  by  Burnikel,  Mehlhorn  and  Schirra  [10].  Fol¬ 
lowing  the  notation  of  Burnikel  [8],  an  incircle  test  is  conveniently  expressed  as  a  quadruple 
(ci ,  02, 03;  04) ,  where  each  Oj  6  {p,  1}  (*  =  !,•••) 4)  is  either  a  point  or  a  line  on  the  plane 
(a  segment  is  seen  by  Burnikel  as  given  by  the  pair  of  its  endpoints  and  by  the  underlying 
line)  and  we  test  whether  04  intersects  the  circle  determined  by  01,03,  and  03. 

The  following  lemma  is  proved  observing  that  the  incircle  test  (pi,P2)P3;P4)  can  be 
answered  by  determining  the  sign  of  4  x  4  determinant  that  is  an  arithmetic  degree  4 
multivariate  polynomial. 
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Lemma  5  [8]  The  incircle  test  (pi,P2,P3;P4)  con  he  solved  with  degree  4. 

Lemma  5  can  be  easily  extended  to  any  dimension  d  >  2.  We  describe  such  test  as 
(pi, . .  where  points  pi, . .  .,pd+i  determine  a  d-dimensional  sphere  and  pd+2 

is  teh  query  point. 

Lemma  6  The  insphere  test  (pi, . .  •  ,Pc<+i  i  Pcf+2)  cny  fixed  dimension  d>2  can  be  solved 

with  degree  d  + 2. 

For  the  construction  of  the  Voronoi  diagram  of  a  set  of  points  and  segments  in  the 
plane  Burnikel  shows  that  the  most  demanding  test  in  terms  of  degree  is  the  incircle  test 

{h,l2J3-,k)  [8]. 

Lemma  7  [8]  The  incircle  test  {li,l2,  h\  h)  con  be  solved  with  degree  40. 

While  the  above  lemmas  provide  an  upper  bound  on  the  degree  of  a  proximity  problem, 
the  next  theorem  gives  a  lower  bound. 

Theorem  1  The  nearest  neighbor  search  problem  for  a  point  set  has  degree  2  in  any  fixed 
dimension  d  >  2. 

Proof:  We  show  the  proof  for  the  case  d  =  2.  The  proof  for  any  other  values  of  d  is 
analogous.  Let  pi  =  P2  =  (3^2)  Pa))  q  =  (xg,yg)  be  three  points  in  the  plane.  In 

order  to  determine  which  of  pi  and  p2  is  the  point  nearest  to  5,  a  point-to-points  distance 
test  must  be  performed. 

This  is  equivalent  to  evaluate  the  sign  of  the  difference  d(pi,g)  -  d(p2,  q)  which,  in  turn, 
is  equivalent  to  the  evaluation  of  the  sign  of  d^{pi,q)  —  d^{p2,q)-  We  claim  that  the  latter 
computation  has  degree  2.  In  fact,  we  show  below  that  d^(pi)9)  ~  d^(P2)9)  cannot  be 
expressed  as  the  product  of  two  degree  1  polynomials,  which  implies  that  the  degree  of  the 
problem  is  2. 

Suppose  that  there  existed  constants  a',  a",  b',  b",  c',  c",  d\  d",  e',  e",  /',  f"  such  that: 

9)  =  a:?  +  Pi  -  2:2  -  P2  -  ^^1^9  +  2x23:,  +  2y2y,  = 

(a'li  -f  b'yi  +  c'x2  +  d'p2  +  +  fyg) '  +  c"x2  +  d"p2  +  e"x,  -f  /"p,)- 

The  above  equality  implies  e'e"  =  0,  since  there  cannot  be  a  term  e'e"x\.  However, 
either  e'  or  e"  is  not  0  because  of  nonzero  terms  having  x,  as  a  factor.  Assume,  w.l.o.g. 
that  e"  ^  0.  Observe  that  d'e"  =  0  because  there  is  no  term  d'e"y2P,;  this  implies  d'  —  0. 
However  we  must  also  have  d'd"  =  -1  because  of  the  term  -p|,  a  contradiction. 

□ 

Observe  that  an  optimal  degree  algorithm  for  the  nearest  neighbor  search  problem  in  a 
planar  point  set  can  be  easily  represented  by  the  brute  force  approach,  where  one  computes 
all  the  distances  between  the  query  point  and  all  other  points  and  reports  the  point  at 
minimum  distance.  However,  such  algorithm  is  both  computationally  inefficient  (it  requires 
quadratic  time)  and  does  not  support  repetitive-mode  queries.  In  Section  3  we  present  an 
optimal  degree  algorithm  whose  query  time  and  space  are  optimal. 
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3  Proximity  Queries  for  Point  Sites  in  the  Plane 

In  this  section,  under  our  standard  assumption  that  all  input  parameters  —  such  as  co¬ 
ordinates  of  sites  and  query  points  —  are  represented  by  h-bit  integers,  we  consider  the 
following  proximity  queries  on  a  set  S  of  point  sites  in  the  plane; 

nearest  neighbor  search:  given  query  point  q,  find  a  site  of  5  whose  Euclidean  distance  from 
q  is  less  than  or  equal  to  that  of  any  other  site; 

k-nearest  neighbors  search:  given  query  point  q,  find  k  sites  of  S  whose  Euclidean  distances 
from  q  are  less  than  or  equal  to  that  of  any  other  site; 

circular  range  search:  given  query  points  q  and  r,  find  the  sites  of  5  that  are  inside  the 
circle  with  center  q  passing  through  r . 

It  is  well  known  that  such  queries  are  efficiently  solved  by  performing  point  location  in 
the  Voronoi  diagram  (possibly  of  higher  order)  V (5)  of  the  sites  [42]. 

We  begin  by  examining  in  Section  3.1  the  geometric  test  primitives  used  by  the  theoreti¬ 
cally  optimal  and  practically  efficient  point  location  methods.  We  identify  three  fundamen¬ 
tal  geometric  test  primitives  for  accessing  the  geometry  of  a  planar  map,  and  we  introduce 
the  concepts  of  '^native”  and  '^ordinary^  point  location  methods. 

In  Section  3.2,  we  show  that  the  “conventional”  approach  of  accessing  the  explicitly  com¬ 
puted  Voronoi  diagram  V(S)  of  the  sites  causes  point-location  queries,  and  hence  proximity 
queries,  to  have  degree  at  least  6. 

In  Sections  3. 3-3.4,  we  describe  our  new  implicit  representation  of  Voronoi  diagrams  for 
point  sites  in  the  plane,  which  allows  us  to  perform  proidmity  queries  with  optimal  degree  2. 

3.1  Test  Primitives  and  Methods  for  Planar  Point  Location 

The  chain  method  [37],  the  bridged  chain  method  [22],  the  trapezoid  method  [41],  the  sub¬ 
division  refinement  method  [35],  and  the  persistent  search  tree  method  [44]  are  popular 
deterministic  techniques  for  point  location  in  a  planar  map  that  combine  theoretical  effi¬ 
ciency  with  good  performance  in  practice  (see,  e.g.,  [21,  42]).  Namely,  denoting  with  n 
the  size  of  the  map,  all  the  above  point  location  methods  have  O(logn)  query  time  and 
0(n  log  n)  preprocessing  time.  The  space  used  is  O(nlogn)  for  the  trapezoid  method  and 
0{n)  for  the  other  methods.  For  monotone  maps,  the  preprocessing  time  is  0{n)  for  the 
chain  method  and  the  bridged  chain  method,  and  O(nlogn)  for  the  other  methods. 

By  a  careful  examination  of  the  query  algorithms  used  in  the  point  location  methods 
presented  in  the  literature,  it  is  possible  to  clearly  separate  the  primitive  operations  that 
access  the  geometry  of  the  map  from  those  that  access  only  the  topology.  We  say  that  a 
point  location  method  is  native  for  a  class  of  maps  if  it  performs  point  locations  queries  in  a 
map  M  of  the  class  by  accessing  the  geometry  of  M  exclusively  through  the  following  three 
geometric  test  primitives  that  discriminate  the  query  point  with  respect  to  the  vertices  and 
edges  of  M : 
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above-below(9,  u)  determine  whether  query  point  q  is  vertically  above  or  below  vertex  v. 

left-right (5,  v)  determine  whether  query  point  q  is  horizontally  to  the  left  or  to  the  right  of 
vertex  v. 

left-right  (5,  e)  determine  whether  query  point  q  is  to  the  left  or  to  the  right  of  edge  e,  this 
operation  assumes  that  edge  e  is  not  horizontal  and  its  vertical  span  includes  q. 

Test  primitive  left-right  (g,  u)  is  typically  used  only  in  degenerate  cases  (e.g.,  in  the  pres¬ 
ence  of  horizontal  edges) . 

A  refinement  of  a  map  M  is  a  map  M'  such  obtained  from  M  by  adding  fictitious 
vertices  and  edges.  Examples  of  refinements  of  a  map  M&TeB,&  triangulation  of  M  and  the 
trapezoidal  decomposition  of  M .  Some  point  location  methods  work  on  refinements  of  the 
original  subdivision  by  means  of  auxiliary  geometric  objects  introduced  in  the  preprocessing 
(e.g.,  triangulation  or  regularization  edges).  We  say  that  a  point  location  method  is  ordinary 
for  a  class  of  maps  if  it  is  native  for  the  refinements  of  the  maps  in  the  class. 

Now,  we  analyze  the  chain  method  [37]  for  point  location  in  a  monotone  map  M.  A 
binary  tree  represents  a  balanced  recursive  decomposition  of  map  M  by  means  of  vertically 
monotone  polygonal  chains  covering  the  edges  of  Af,  called  separators.  A  point  location 
query  consists  of  traversing  a  root-to-leaf  path  in  this  tree,  where  at  each  node  we  determine 
whether  the  query  point  q  is  to  the  left  or  to  right  of  the  separator  associated  with  the  node. 
The  discrimination  of  point  q  with  respect  to  a  separator  a  is  performed  in  two  steps; 

1.  we  find  the  edge  e  of  cr  whose  vertical  span  includes  point  q  by  means  of  binary  search 
on  the  y  coordinates  of  the  vertices  of  <t,  which  consists  of  performing  a  sequence  of 
a  logarithmic  number  of  above-below(g,  u)  tests; 

2.  we  discriminate  q  with  respect  to  a  by  performing  test  left-right(g,  e). 

In  the  special  case  that  separator  a  has  horizontal  edges,  the  discrimination  of  point  q 
with  respect  to  a  uses  also  test  primitive  left-right (g,  u).  Hence,  the  chain  method  is  native 
for  monotone  maps.  For  a  map  M  that  is  not  monotone,  fictitious  “regularization  edges  are 
added  to  M  and  point  location  in  M  is  reduced  to  point  location  in  the  resulting  refinement 
M'  of  M.  Hence,  the  chain  method  is  ordinary  for  general  maps. 

In  the  bridged  chain  method  [22],  the  technique  of  fractional  cascading  [15, 16]  is  applied 
to  the  sets  of  y-coordinates  of  the  separators.  Fractional  cascading  establishes  “bridges” 
between  the  separator  of  a  node  and  the  separators  of  its  children  such  that  there  are 
0(1)  vertices  between  any  two  consecutive  bridges.  Hence,  except  for  the  separator  of  the 
root.  Step  1  can  be  executed  with  0(1)  above-below(g,  v)  tests  for  the  vertices  between  two 
consecutive  bridges.  The  bridged  chain  method  is  ordinary  for  general  maps  and  native  for 
monotone  maps. 

A  similar  analysis  shows  that  all  efficient  point  location  methods  described  in  the  liter¬ 
ature  are  ordinary  for  general  maps.  More  specifically,  we  have; 
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Lemma  8  The  trapezoid  method  and  the  persistent  search  tree  method  are  native  for  general 
maps.  The  chain  method  and  the  bridged  chain  method  are  ordinary  for  general  maps  and 
native  for  monotone  maps.  The  subdivision  refinement  method  is  ordinary  for  general  maps 
and  is  native  for  trivial  maps  consisting  of  a  single  triangle. 

Hence,  all  the  known  planar  point  location  methods  described  in  the  literature  are 
ordinary  for  Voronoi  diagrams. 

3.2  Explicit  Voronoi  Diagrams 

Let  5  be  a  set  of  n  point  sites  in  the  plane,  where  each  site  is  a  primitive  point  with  6-bit 
integer  coordinates.  The  Voronoi  diagram  y(5)  of  S  is  said  to  be  explicit  if  the  coordinates 
of  the  vertices  of  V{S)  are  computed  and  stored  with  exact  arithmetic,  i.e.,  as  rational 
numbers  (pairs  of  integers). 

Lemma  9  The  left-right (g,  e)  test  primitive  in  an  explicit  Voronoi  diagram  of  point  sites 
in  the  plane  has  degree  6. 

Proof:  Let  e  =  (ui,U2)  be  a  Voronoi  edge  such  that  ui  =  (a:i,  j/i)  is  equidistant  from  three 
sites  a  =  {xa,ya),  b  =  (X6,yb),  c=  (x^yc)  and  uj  =  {x2,y2)  is  equidistant  from  three  sites 
6  =  {xb^Vb),  c  =  (xc.yc),  and  d  =  {xd,yd)-  See  Figure  1.  In  an  explicit  Voronoi  diagram, 
te.st  primitive  left-right  (g,  e)  that  determines  whether  query  point  g  =  (xq^yg)  is  to  the  left 
or  to  the  right  of  edge  e  =  (‘^11^2)  is  equivalent  to  evaluating  the  sign  of  the  following 
determinant: 


Xq 

Vq 

9 

Vq 

1 

Xq 

yg 

1 

A  = 

Xi 

yi 

1 

JVi 

Y\ 

2Wi 

1 

1 

”  2W1W2 

Xl 

Yi 

Wi 

X2 

y2 

1 

2V^ 

2W2 

1 

X2 

Y2 

W2 

where 


xl  +  yl  Pa  1 

Xa  xI  +  pI  1 

Xa  Pa  1 

xl  +  yb  Pb  1 

,  Vi  = 

Xb  xl  +  yl  1 

,  Wi  = 

Xb  Pb  1 

xl  +  yl  Pc  1 

Xc  xl  +  yl  1 

Xc  Pc  1 

and  X2,Y2,  and  W2  have  similar  expressions  obtained  replacing  in  the  above  determinants 
Xc  with  Xd  and  Pc  with  yj.  Evaluating  the  sign  of  A  is  equivalent  to  evaluating  the  signs  of 
Wu  W2  and  of  A'. 

By  using  the  notation  introduced  in  Section  2.1,  we  can  rewrite  Xi  and  Vj  as  a^,  and  W, 
as  (t  =  1,2).  Hence  A'  is  a  degree-6  multivariate  polynomial  since  it  can  be  rewritten 
as 

-  a^cr^)  -  Q;(cr^(v^  -  a^a^)  +  — ^(2.3,4)  q,6  ^  q,6  — ^(3)  ^,6 
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Figure  1;  Illustration  for  Lemma  9. 


An  algorithm  for  proximity  queries  on  a  set  S  of  point  sites  in  the  plane  is  said  to 
be  conventional  if  it  computes  the  explicit  Voronoi  diagram  V (5)  of  5  and  then  performs 
point  location  queries  on  V  (S)  with  an  ordinary  method.  Note  that  the  class  of  conventional 
proximity  query  algorithms  includes  all  the  efficient  algorithms  presented  in  the  literature. 
A  conventional  proximity  query  algorithm  needs  to  perform  test  primitive  left-right(g,  e). 
Thus,  by  Lemma  9  we  have: 

Theorem  2  Conventional  algorithms  for  the  follovoing  proximity  query  problems  on  a  set 
of  point  sites  in  the  plane  have  degree  at  least  6: 

•  nearest  neighbor  query; 

•  k -nearest  neighbor  query; 

•  circular  range  query. 

We  observe  that  a  degree-6  algorithm  implies  that  a  fc-bit  arithmetic  unit  can  handle 
with  native  precision  queries  for  points  in  a  grid  of  size  at  most  2^^®  X  2^^^.  For  example, 
if  =  32,  the  points  that  can  be  treated  with  single-precision  arithmetic  belong  to  a  grid 
of  size  at  most  64  x  64. 

3.3  Implicit  Voronoi  Diagrams 

Let  S  be  a  set  of  n  point  sites  in  the  plane,  and  recall  our  assumption  that  each  site  or 
query  point  is  a  primitive  point  with  6-bit  integer  coordinates.  We  say  that  a  number  s 
is  a  semi-integer  if  it  is  a  rational  number  of  the  type  s  =  m/2  for  some  integer  m.  The 
implicit  Voronoi  diagram  V*(5)  of  5  is  a  representation  of  the  Voronoi  diagram  V(5)  of  5 
that  consists  of  a  topological  component  and  of  a  geometric  component.  The  topological 
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component  of  V*{S)  is  the  planar  embedding  of  1^(5),  represented  by  a  suitable  data  struc¬ 
ture  (e.g.,  doubly-connected  edge  lists  [42]  or  the  data  structure  of  [30]).  The  geometric 
component  of  y*(5)  stores  the  following  geometric  information  for  each  vertex  and  edge  of 
the  embedding: 


•  For  each  vertex  v  of  V(5),  V*(5)  stores  semi-integers  x*(u)  and  y*{v)  that  approx¬ 
imate  the  X-  and  ^-coordinates  y{v)  of  v,  We  provide  the  definition  of  y*{v)  below. 
The  definition  of  x*{v)  is  analogous. 


y(u) 

L2/(«)J  +  2 

ot>  _  1 
^  2 
0 


0  <  y(u)  <  2*’  -  1,  y{v)  integer 
0  <  3/{u)  <  2*’  -  1,  y{v)  not  integer 
y{v)  >  2*’  -  1 
y(u)  <  0 


Note  that  semi-integers  x*(t;)  and  y*(t;)  can  be  stored  with  (6-f  l)-bits. 

•  For  each  non-horizontal  edge  e  of  V(5),  V'(5)  stores  the  pair  of  sites  e{e)  and  r(e) 
such  that  e  is  a  portion  of  the  perpendicular  bisector  of  £{e)  and  r(e),  and  £{e)  is  to 
the  left  of  r(e). 


Equipped  with  the  above  information,  the  three  test  primitives  for  point  location  can 
be  performed  in  the  implicit  Voronoi  diagram  V*(5)  as  follows: 


above-below(y,  v)  compare  the  y-coordinate  of  q  with  y*(n); 
left-right  (g,u)  compare  the  x-coordinate  of  q  with  x*(u); 

left-right(g,  e)  compare  the  Euclidean  distances  of  point  q  from  sites  ^(e)  and  r(c). 

Since  the  query  point  q  is  by  assumption  a  primitive  point  whose  coordinates  are  &-bit 
integers,  we  have  that  y{q)  <  y{v)  if  and  only  if  y(g)  <  y*{v),  where  testing  the  latter 
inequality  has  degree  1.  Similar  considerations  apply  to  testing  x{q)  <  x{v).  This  proves 
the  correctness  of  our  implementation  of  above-below(g,t;)  and  left-right(g,  u). 

The  correctness  of  the  above  implementation  of  test  left-right  (g,  e)  follows  directly  from 
the  definition  of  Voronoi  edges.  Thus,  in  an  implicit  Voronoi  diagram,  test  left-right(g,e) 
can  be  implemented  with  a  point-to-points  distance  test  that  has  degree  2  (Lemma  4). 

Hence,  we  obtain  the  following  lemmas: 


Lemma  10  Test  primitives  above-below(g,  v)  and  left-right(g,  u)  in  an  implicit  Voronoi  di¬ 
agram  of  point  sites  in  the  plane  can  be  performed  in  0(1)  time  and  with  degree  1. 

Lemma  11  Test  primitive  left-right  (g,e)  in  an  implicit  Voronoi  diagram  of  point  sites  in 
the  plane  can  be  performed  in  0(1)  time  and  with  degree  2. 
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In  order  to  execute  a  native  point  location  algorithm  in  an  implicit  Voronoi  diagram, 
we  only  need  to  redefine  the  implementation  of  the  three  test  primitives.  By  having  encap¬ 
sulated  the  geometry  in  the  test  primitives,  no  further  modifications  are  needed.  Hence,  by 
Lemmas  10-11  we  obtain: 

Lemma  12  For  any  native  method  on  a  class  of  maps  that  includes  Voronoi  diagrams,  a 
point  location  query  in  an  implicit  Voronoi  diagram  has  optimal  degree  2  and  has  the  same 
asymptotic  time  complexity  as  a  point  location  query  in  the  corresponding  explicit  Voronoi 
diagram. 

In  order  to  compute  the  implicit  Voronoi  diagram  V’*(5),  we  begin  by  constructing  the 
Delaunay  triangulation  of  S,  denoted  DT{S),  by  means  of  the  O  (n  log  n)-time  algorithm 
of  [30],  which  has  degree  4  because  its  most  expensive  operation  in  terms  of  the  degree  is  the 
incircle  test  (see  Lemma  5).  The  topological  structure  of  V(iS)  and  the  sites  ^(c)  and  r(e) 
for  each  edge  e  of  V{5)  are  immediately  derived  from  DT{S)  by  duality.  Next,  we  compute 
the  approximations  x‘(u)  and  y*{v)  for  each  vertex  v  oiV{S)  by  means  of  integer  division. 
Let  a,  6,  and  c  be  the  three  sites  of  S  that  define  vertex  v.  Adopting  the  same  notation  as 
in  the  proof  of  Lemma  9,  the  y-coordinate  y{v)  of  v  is  given  by  the  ratio  y{v)  —  where 
Yi  is  a  variable  of  arithmetic  degree  3  and  W2  is  a  variable  of  arithmetic  degree  2,  and 
similarly  for  x(u).  Hence,  the  computation  of  and  y*{v)  has  degree  3.  We  summarize 
the  above  analysis  as  follows. 

Lemma  13  The  implicit  Voronoi  diagram  of  n  point  sites  in  the  plane  can  be  computed  in 
0(n  log n)  time,  0{n)  space,  and  with  degree  4- 

Theorem  3  Let  S  be  a  set  of  n  point  sites  in  the  plane.  There  exists  an  0(n) -space  data 
structure  for  S,  based  on  the  implicit  Voronoi  diagram  V*{S),  that  can  be  computed  in 
O(nlogn)  time  with  degree  5,  and  supports  nearest  neighbor  queries  in  O(logn)  time  with 
optimal  degree  2. 

Proof:  We  perform  point  location  in  the  implicit  Voronoi  V*(5)  diagram  of  5  using  a 
native  method  for  monotone  maps  with  optimal  space  and  query  time,  such  as  the  bridged- 
chain  method  or  the  persistent  search  tree  method.  The  space  requirement  and  the  query 
degree  and  time  follow  from  the  performance  of  these  methods  and  from  Lemma  12. 

Regarding  the  preprocessing  time,  by  Lemma  13,  the  construction  of  the  implicit  Voronoi 
V(S)  takes  time  O(nlogn)  time  with  degree  4.  In  order  to  construct  the  point  location 
data  structure,  we  also  need  an  additional  test  primitive  that  consists  of  comparing  the 
y-coordinates  of  two  Voronoi  vertices.  E.g.,  this  primitive  is  used  to  establish  bridges  in 
the  bridged-chain  method  (see  Section  3.1)  and  to  sort  the  vertices  by  y-coordinate  in  the 
persistent  location  method.  We  can  show  that  this  primitive  has  degree  5.  □ 
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3.4  Implicit  Higher  Order  Voronoi  Diagrams 

In  this  section,  we  introduce  implicit  higher  order  Voronoi  diagrams  for  point  sites  in  the 
plane,  and  extend  the  results  of  Section  3.3  to  fc-nearest  neighbors  and  circular  range  search 
queries. 

The  definition  of  the  implicit  order-k  Voronoi  diagram  (S)  of  set  5  of  point  sites  in 
the  plane  is  analogous  to  that  given  in  Section  3.3  for  Voronoi  diagrams.  A  vertex  v  of  Vi- (5) 
is  represented  by  its  approximate  coordinates  a:*(v)  and  y*{v),  and  a  non-horizontal  edge  e 
of  Vk{S)  stores  the  pair  of  sites  £{e)  and  r(e)  such  that  e  is  a  portion  of  the  perpendicular 
bisector  of  ({e)  and  r(e),  and  £{e)  is  to  the  left  of  r(e). 

Lemmas  10-11  immediately  hold  also  for  Vi;(S),  and  we  obtain: 

Lemma  14  For  any  native  method  for  monotone  maps,  a  point  location  Query  in  an  implicit 
order-k  Voronoi  diagram  has  optimal  degree  2  and  has  the  same  asymptotic  time  complexity 
as  a  point  location  query  in  an  explicit  order-k  Voronoi  diagram. 

The  order-/:  Voronoi  diagram  Vk{S)  for  a  set  S  of  n  point  sites  has  0{k{n  —  k))  vertices, 
edges,  and  faces,  and  can  be  obtained  from  the  order  fc  —  1  implicit  Voronoi  diagram  T4_i  (5) 
by  intersecting  each  face  of  V^—i  (S)  with  the  (order  1)  Voronoi  diagram  of  a  suitable  subset 
of  the  vertices  of  5  [36].  As  shown  in  [36, 14],  V*:(S)  can  be  computed  in  O (k (n-k)y/n  log  n) 
time.  Since  the  construction  is  based  on  iteratively  computing  Voronoi  diagrams  by  using 
the  incircle  test,  which  is  the  most  expensive  operation  in  terms  of  degree,  the  overall  degree 
of  the  preprocessing  is  4  (Lemma  5).  Hence,  we  obtain. 

Lemma  15  The  implicit  order-k  Voronoi  diagram  of  n  point  sites  in  the  plane  can  be 
computed  in  0{k{n-  k)y/n\og  n)  time,  0{k{n-k))  space,  and  with  degree  4. 

Point  location  in  the  order-A:  Voronoi  diagram  solves  A:-nearest  neighbors  queries.  Hence, 
by  Theorem  1  and  Lemmas  14-15,  we  obtain: 

Theorem  4  Let  S  be  a  set  of  n  point  sites  in  the  plane  and  k  an  integer  with  1  <  k  <  n-1. 
There  exists  an  0{k{n-k))-space  data  structure  for  S,  based  on  the  implicit  order-k  Voronoi 
diagram  Vi^{S),  that  can  be  computed  inO{k{n-k)y/nlogn)  time  with  degree  b  and  supports 
k-nearest  neighbors  queries  in  0(logn-l-  k)  time  with  optimal  degree  2. 

Circular  range  search  queries  in  a  set  5  of  n  point  sites  can  be  reduced  to  a  sequence 
of  2’-nearest  neighbors  queries  in  V2'(‘5)>  »  =  0,  •  •  •,logn  [7].  This  approach  yields  a  data 
structure  with  O(n^)  space  and  preprocessing  time,  and  O (log  n  log  log  n  4-  A:)  query  time, 
where  k  is  the  size  of  the  output.  Hence,  with  analogous  reasoning  as  above,  we  obtain  the 
following  theorem. 

Theorem  5  Let  S  be  a  set  of  n  point  sites  in  the  plane.  There  exists  an  0{n^)-space  data 
structure  for  S,  based  on  implicit  order-k  Voronoi  diagrams,  that  can  be  computed  in  0(n^) 
time  with  degree  5  and  supports  circular  range  search  queries  in  O  (log  n  log  log  n  -1-  A:)  time 
with  optimal  degree  2. 
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The  space  and  preprocessing  time  of  Theorems  4-5  and  the  query  time  of  Theorem  - 
5  can  be  substantially  improved  while  preserving  the  same  degree  bounds  using  the  data 
structures  presented  in  [1,  2,  13]. 


4  Proximity  Queries  for  Point  Sites  in  3D  Space 

In  this  section,  we  consider  the  following  proximity  query  on  a  set  5  of  point  sites  in  three- 
dimensional  (3D)  space: 

nearest  neighbor  search:  given  query  point  q,  find  a  site  of  5  whose  Euclidean  distance  from 
q  is  less  than  or  equal  to  that  of  any  other  site; 

We  recall  our  assumption  that  the  sites  and  query  points  are  primitive  points  represented 

by  6-bit  integers.  r  •  •  + 

As  for  the  two-dimensional  case,  such  query  is  efficiently  answered  by  performing  pomt 

location  in  the  3D  Voronoi  diagram  of  S.  Test  primitives  and  methods  for  spatial  point 
location  are  described  in  Section  4.1.  Section  4.2  shows  that  “conventional  algorithms 
require  degree  8.  A  degree  3  algorithm  based  on  “implicit”  3D  Voronoi  diagrams  is  then 

given  in  Section  4.3. 

4.1  Test  Primitives  and  Methods  for  Spatial  Point  Location 

There  are  only  two  known  efficient  spatial  point  location  methods  for  cell-complexes  that  are 
applicable  to'3D  Voronoi  diagrams:  the  separating  surfaces  method  [12,  47],  which  extends 
the  chain-method  [37],  and  the  persistent  planar  location  method  [43],  which  extends  the 
persistent  search  tree  method  [44].  Let  N  be  the  number  of  facets  of  a  cell-complex  C. 
The  query  time  is  0{\og^  N)  for  both  methods.  The  space  used  is  0{N)  for  the  separating 
surfaces  method  and  OCVlog^  N)  for  the  persistent  planar  location  method.  Both  methods 
are  restricted  to  convex  cell-complexes.  The  separating  surffices  method  is  further  restricted 
to  acyclic  convex  cell-complexes,  where  the  dominance  relation  among  cells  in  the  2-direction 

is  scvclic. 

As  in  Section  3.1,  we  can  separate  the  primitive  operations  that  access  the  geometry  of 
the  cell-complex  from  those  that  access  only  the  topology.  We  say  that  a  point  location 
method  is  native  for  a  class  of  3D  cell-complexes  if  it  performs  point  locations  queries  in  a 
cell  complex  C  of  the  class  by  accessing  the  geometry  of  C  exclusively  through  the  following 
three  geometric  test  primitives  that  discriminate  the  query  point  with  respect  to  the  vertices 

and  edges  of  C: 

above-below(g,u)  compare  the  z-coordinate  of  the  query  point  q  with  the  z-coordinate  of 
vertex  u. 

left-right  (g,  v)  compare  the  i-coordinate  of  the  query  point  q  with  the  i-coordinate  of  ver¬ 
tex  V. 
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front-rear(g,u)  compare  the  t/-coordinate  of  the  query  point  q  with  the  y-coordinate  of 
vertex  v. 

left-right  e^y)  compare  the  ly-projection  Qxy  of  the  query  point  q  with  the  xy-projection 
of  edge  txy  This  operation  assumes  that  e^y  is  not  parallel  to  the  x-axis  and  its  y-span 
includes  qxy. 

above-below(gyj,eyj)  compare  the  yx-projection  qyz  of  the  query  point  q  with  the  yz- 
projection  of  edge  Cy*.  This  operation  assumes  that  eyz  is  not  parallel  to  the  y-axis 
and  its  z-span  a  includes  qyz . 

above-below(g,  /)  determine  whether  query  point  q  is  above  or  below  a  facet  /;  this  opera¬ 
tion  assumes  that  facet  /  is  not  parallel  to  the  x-axis  and  that  the  xy-projection  of  / 
contains  the  xy-projection  of  q. 

Test  primitives  front-rear(g,  v)  and  left-right(g,  v)  are  used  only  in  degenerate  cases  (e.g. 
in  the  presence  of  facets  parallel  to  the  z-axis  and  in  cases  where  Cxy  is  horizontal). 

Now,  we  analyze  the  separating  surfaces  method  for  spatial  point  location  [12,  47]  in 
acyclic  cell-complexes.  Separating  surfaces  are  the  3D  analogue  of  separators  of  monotone 
maps.  Their  existence  is  guaranteed  by  the  acyclicity  of  the  cell-complex.  Thus,  a  point 
location  query  consists  of  traversing  a  root-to-leaf  path  in  the  separating  surface  tree,  where 
at  each  node  we  determine  whether  the  query  point  q  is  to  above  or  below  the  separating 
surface  associated  with  the  node.  The  discrimination  of  point  q  with  respect  to  a  separating 
a  is  performed  in  two  steps: 

1.  by  means  of  a  planar  point  location  query  for  the  xy-projection  qxy  of  q  in  the  xy  pro¬ 
jection  of  <7,  we  find  the  facet  /  of  (T  whose  xy  projection  contains  qxy  If  an  ordinary 
point  location  method  is  used,  this  step  uses  primitives  front-rear(g,  v),  left-right (g,  v), 
and  left-right(gxy,exy). 

2.  we  discriminate  g  with  respect  to  tr  by  performing  test  above-below(g,  /). 

In  the  special  cases  that  cell-complex  C  has  facets  parallel  to  the  z-axis,  the  discrim¬ 
ination  of  point  g  with  respect  to  cr  uses  also  test  primitives  above-below(g,  u).  Thus,  the 
separating  surfaces  method  is  native  for  acyclic  convex  cell-complexes. 

A  similar  analysis  shows  that  also  the  persistent  planar  location  method  is  native  for 
convex  cell-complexes.  More  specifically,  we  have: 

Lemma  16  The  separating  surfaces  method  is  native  for  acyclic  convex  cell-complexes. 
The  persistent  planar  location  method  is  native  for  convex  cell-complexes. 

Hence,  all  the  known  spatial  point  location  methods  described  in  the  literature  are  native 
for  3D  Voronoi  diagrams. 
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4.2  Explicit  Voronoi  Diagrams 

Let  5  be  a  set  of  n  point  sites  in  3D,  where  each  site  is  a  primitive  point  with  h-bit  integer 
coordinates.  The  3D  Voronoi  diagram  V{S)  of  5  is  said  to  be  explicit  ii  the  coordinates 
of  the  vertices  of  y(S)  are  computed  and  stored  with  exact  arithmetic,  i.e.,  as  rational 
numbers  (pairs  of  integers). 

Lemma  17  The  left-right  (g^y,  e^y)  test  primitive  in  an  explicit  Voronoi  diagram  of  point 
sites  in  SD  space  has  degree  8. 

Proof:  The  reasoning  is  the  same  as  in  the  proof  of  Lemma  9.  Let  ex,y  =  (^^11^2))  where  ui 
and  V2  are  the  ly-projections  of  two  adjacent  vertices  u  and  v  of  V(5);  let  u^e  equidistant 

from  the  four  primitive  sites  a  =  {Xa^Va)^  ^  =  {xb^Vb)^  ^  d  =  an 

V  from  «  =  (X..!,.),  b  =  c  =  (x„!,,).  and  h  =  (xn,!/n).  A  conventional  algonthm 

determines  whether  query  point  ?  =  (i,,!;,)  is  to  the  left  or  to  the  right  of  the  oriented 
edge  e  =  (ui ,  V2)  by  evaluating  the  sign  of  the  following  determinant: 

Xq  Vg  g  Xq  Pg  I  Xg  Pg  1 

A  1  —  1  =  ~  -JU',  u'o  1 

Xi  Y2  W2 


where 


A'i  = 


Xg  Pg 

q 

A  = 

Xi  Pi 

1 

zzi 

X2  P2 

1 

x^-^r  p’i'V  4 

ya 

xl  +  Pb  +  4 

yb 

Zb 

xl+Pc+  4 

yc 

Zc 

Xd  +  yd  +  4 

yd 

Zd 

Xg  Pg 

X,  Y, 

2W'l  2My 

^  2^ 


A' 

2W1W2  ’ 


xI  +  pI  +  4 
xl  +  yl  +  4 

X^  +  Pc  +  4 

Xd  +  yd  +  4 


Wi  = 


O-JLiVl  4^  **  Z  ^  i * 

Xd  with  Xh,  yd  with  ph,  and  Zd  with  zh-  _  r  Tirr  iirr  j  f  a/ 

Evaluating  the  sign  of  A  is  equivalent  to  evaluating  the  signs  of  Wi,  W2  and  01  Is  .  ^ 

By  using  the  notation  introduced  in  Section  2.1,  we  can  rewrite  X,-  and  Yi  as  a  ,  and  Wi 
q3  _  12).  Hence  A'  is  a  degree-8  multivariate  polynomial  since  it  can  be  rewritten 

as 

-  a^'ci^)  -  a(aW  -  a'^a®) -b  a'‘a3  -  aW  —>(2,3,4)  — >(^)  a®. 


An  algorithm  for  nearest  neighbor  queries  on  a  set  5  of  point  sites  in  3D  space  is 
said  to  be  conventional  if  it  computes  the  explicit  3D  Voronoi  diagram  V(S)  of  S  and 
then  performs  point  location  queries  on  V(5)  with  a  native  method.  Recall  that  the  class 
of  conventional  nearest  neighbor  query  algorithms  includes  the  two  efficient  algorithms 
presented  in  the  literature.  A  conventional  proximity  query  algorithm  needs  to  perform 
test  primitive  left-right Thus,  by  Lemma  17,  we  have; 

Theorem  6  Conventional  algorithms  for  the  nearest  neighbor  querp  problem  on  a  set  of 
point  sites  in  3D  space  have  degree  at  least  8. 
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4.3  Implicit  Voronoi  Diagrams 

The  definition  of  the  implicit  3D  Voronoi  diagram  V*(5)  of  a  set  of  5  of  point  sites  in  3D 
space  is  a  straightforward  extension  of  the  definition  for  two-dimensional  Voronoi  diagrams 
given  in  Section  3.3.  Namely  V*{S)  stores  the  topological  structure  of  the  3D  Voronoi 
diagram  V (5)  of  5  (e.g.,  the  data  structure  of  [20])  and  the  following  geometric  information 
for  each  vertex  and  facet: 


•  For  each  vertex  v  of  V{S),  V*{S)  stores  the  semi-integer  (6-t-  l)-bit  approximations 
x'‘{v),  y*{v)  and  2*(u)  of  the  x-,  y-,  and  z-coordinates  of  v. 

•  For  each  facet  /  of  V’(S)  that  is  not  parallel  to  any  of  three  Cartesian  planes,  V*(5) 
stores  the  pair  of  sites  £{f)  and  r(/)  such  that  /  is  a  portion  of  the  perpendicular 
bisector  of  £{f)  and  r(/),  and  £{f)  is  below  r(/). 

The  tests  above-below(g,  u),  left-right(5,u),  front-rear(?,  u)  can  be  implemented  compar¬ 
ing  the  X-,  y-  and  z-coordinate  of  query  point  q  with  y(u)*,  and  x(u)*  respectively. 

With  the  same  reasoning  as  for  the  two-dimensional  case  (see  Section  3.3),  it  is  easy  to  see 
that  such  implementations  are  correct. 

Lemma  18  Test  primitives  above-below(g,  u),  left-right(g,  v),  front-rear(g,  u)  in  an  implicit 
Voronoi  diagram  of  3D  point  sites  can  be  performed  in  0(1)  time  and  with  degree  1. 

Test  primitive  above-below(g, /)  is  implemented  by  comparing  the  Euclidean  distances 
of  point  q  from  the  two  sites  £{e)  and  r(e)  of  which  /  is  the  perpendicular  bisector,  with 
a  point-to-points  distance  test.  The  implementation  is  correct  by  the  definition  of  Voronoi 
facet.  Thus,  by  Lemma  4,  we  have. 

Lemma  19  Test  primitive  above-below(9, /)  in  an  implicit  Voronoi  diagram  of  SD  point 
sites  can  be  performed  in  0{1)  time  and  with  degree  2. 


Finally,  test  left-right(ga:y,  ea,y)  is  implemented  by  determining  the  sign  of  the  equation 
of  the  line  that  contains  edge  e^y  when  computed  at  point  g^y. 

Lemma  20  Test  primitive  left-right (g^y,  e^y)  in  an  implicit  Voronoi  diagram  of  SD  point 
sites  can  be  performed  in  0{l)  time  and  with  degree  3. 


Proof:  The  line  containing  the  oriented  edge  e^y  is  the  projection  on  the  xy-plane  of  the 
intersection  of  two  planes  containing  two  facets  of  the  three  dimensional  Voronoi  diagram. 
Let  a,x  -t-  biy  +  CiZ  +di  =  0  be  the  equation  of  one  such  plane  (i  =  1, 2).  The  projection  of 
their  intersection  on  the  ly-plane  is 
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02 


Cl 


X  + 
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di  Cl 
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Test  left-right  (gi-y,  e^y)  consists  of  determining  the  sign  of 


Oi  Cl 

02  C2 


Xa  + 


bi  Cl 

h  C2 


Vg  + 


di  Cl 
d2  C2 


which  is  a  multivariate  polynomial  having  arithmetic  degree  3,  since  it  can  be  rewritten  as 

aa^  -h  aa^  + 


k.2  I  ^,^,2  1  ^,3  _ 


In  order  to  execute  a  native  point  location  algorithm  in  an  implicit  3D  Voronoi  diagram, 
we  only  need  to  redefine  the  implementation  of  the  five  test  primitives.  By  having  encap¬ 
sulated  the  geometry  in  the  test  primitives,  no  further  modifications  are  needed.  Hence,  by 
Lemmas  18-20  we  obtain; 

Lemma  21  For  any  native  method  on  a  class  of  cell-complexes  that  includes  SD  Voronoi 
diagrams,  a  point  location  query  in  an  implicit  3D  Voronoi  diagram  has  degree  S  and  has 
the  same  asymptotic  time  complexity  as  a  point  location  query  in  an  explicit  3D  Voronoi 
diagram. 

The  Voronoi  diagram  of  n  point  sites  in  3D  space  is  an  acyclic  convex  cell  complex  with 
N  =  0(n^)  facets.  Hence,  using  the  separating  surfaces  method  on  the  implicit  3D  Voronoi 
diagram  yields  the  following  result: 

The  implicit  Voronoi  diagram  V*(5)  of  a  set  5  of  n  points  in  3D  space  can  be  constructed 
by  computing  the  3D  Delaunay  triangulation  with  the  incremental  algorithm  by  Joe  [33], 
whose  time  complexity  and  storage  is  0{n^)  (see  also  [40]).  Since  the  most  demanding 
operation  of  the  algorithm  in  terms  of  degree  is  the  3D  insphere  test,  from  Lemma  6  we 
have  that  the  degree  of  the  cdgorithm  that  computes  V (S)  is  5.  As  in  the  planar  case, 
the  topological  structure  of  V(S')  and  the  sites  £(f)  and  ?*(/)  for  each  edge  e  of  V(S')  are 
immediately  derived  from  DT(S)  by  duality.  Omitting  various  details,  we  have: 

Lemma  22  The  implicit  Voronoi  diagram  of  a  set  of  n  point  sites  in  3D  space  can  be 
computed  in  O(n^)  time  and  space,  and  with  degree  5. 

Lemmas  21  and  22  lead  to  the  following  theorem. 

Theorem  7  Let  S  be  a  set  of  n  point  sites  in  3D  space.  There  exists  an  0{n^)-space  data 
structure  for  S  that  can  be  computed  in  O(n^)  time  with  degree  7  and  supports  nearest 
neighbor  queries  in  O  (log^  n)  time  with  degree  3 . 

Note  that  the  degree  7  of  the  preprocessing  is  due  to  an  additional  test  primitive  that 
consists  of  comparing  the  y-coordinates  of  two  Voronoi  vertices 

Although  the  algorithm  for  nearest  neighbor  queries  proposed  in  this  section  has  nonop- 
timal  degree  3,  it  is  a  practical  approach  for  the  important  application  scenario  where  the 
primitive  points  are  pixels  on  a  computer  screen.  On  a  typical  screen  with  about  2^®  x  2^® 
pixels,  our  nearest  neighbor  query  can  be  executed  with  the  standard  integer  arithmetic  of 
a  32-bit  processor. 
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5  Proximity  Queries  for  Point  and  Segment  Sites  in  the 
Plane 

In  this  section,  we  consider  the  following  proximity  query  on  a  set  5  of  point  and  segment 
sites  in  the  plane: 

nearest  neighbor  search:  given  query  point  g,  find  a  site  of  5  whose  Euclidean  distance  from 
q  is  less  than  or  equal  to  that  of  any  other  site. 

As  for  the  other  queries  studied  in  the  previous  sections,  such  query  is  efficiently  solved 
by  performing  point  location  in  the  Voronoi  diagram  of  the  set  of  point  and  segment 
sites  [42]. 

The  test  primitives  needed  by  such  approach  are  described  in  Section  5.1.  Section  5.2 
shows  that  the  “conventional”  approach  requires  degree  64.  A  degree  6  algorithm  based  on 
“implicit”  Voronoi  diagrams  is  then  given  in  Section  5.3. 

5.1  Test  Primitives  and  Methods 

The  Voronoi  diagram  V (5)  of  a  set  5  of  point  and  segment  sites  is  a  map  whose  edges  are 
either  straight-line  segments  or  arcs  of  parabolas.  Hence,  in  general  V (5)  is  neither  convex 
nor  monotone.  In  order  to  perform  point  location  in  V(S),  we  refine  V{S)  into  a  map 
with  monotone  edges  as  follows.  If  edge  e  of  V(S)  is  an  arc  of  parabola  whose  point  p  of 
maximum  (or  minimum)  y-coordinate  is  not  a  vertex,  we  split  e  into  two  edges  by  inserting 
a  fictitious  vertex  at  point  p.  We  call  the  resulting  map  the  extended  Voronoi  diagram  V'{S) 
of  5. 

The  persistent  search  tree  method  and  the  trapezoid  method  can  be  used  as  native 
methods  on  the  extended  Voronoi  diagram,  where  the  test  primitives  are  the  same  as  those 
defined  in  Section  3.1  for  point  sites.  If  we  want  to  use  the  chain  method  or  the  bridged 
chain  method,  we  need  to  do  a  further  refinement  that  transforms  the  map  into  a  monotone 
map  by  adding  vertical  fictitious  edges  emanating  from  the  fictitious  vertices  previously 
inserted  along  the  parabolic  edges. 

Lemma  23  The  trapezoid  method  and  the  persistent  search  tree  method  are  native,  and  the 
chain  method  and  the  bridged  chain  method  are  ordinary  for  extended  Voronoi  diagrams  of 
point  and  segment  sites. 

5.2  Explicit  Voronoi  Diagrams 

Let  S  be  a  set  of  n  points  and  segment  sites.  The  extended  Voronoi  diagram  V’'(5)  of  S  is 
said  to  be  explicit  if  the  coordinates  of  the  vertices  of  V'(S)  are  computed  and  stored  with 
exact  arithmetic,  i.e.,  as  algebraic  numbers  [9,  49]. 

In  the  following  lemma,  we  analyze  the  degree  of  test  primitive  left-right  (g,e)  for  a 
straight-line  edge  e  of  an  explicit  extended  Voronoi  diagram. 
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Lemma  24  The  left-right (g,  e)  test  primitive  for  a  straight-line  edge  e  in  an  explicit  ex¬ 
tended  Voronoi  diagram  of  point  and  segment  sites  in  the  plane  has  degree  64. 

Proof:  Let  e  =  (ui,V2),  such  that  vi  =  (xi,yi)  is  equidistant  from  three  segments  si,  S2, 
and  S3  and  V2  is  from  three  segments  si,  S2,  and  S4.  See  Figure  2. 


Figure  2:  Illustration  for  Lemma  24. 


We  show  that  the  test  left-right(5,  e)  for  determining  whether  the  query  point  q  =  (xg,  y,) 
is  to  the  left  or  to  the  right  of  (ui,V2)  has  degree  64.  Namely,  let  +  biy  +  Ci  =  0 
(?.  =  1,2, 3, 4)  be  the  equation  of  the  line  containing  segment  s^.  A  conventional  algorithm 
computes  the  above  test  by  evaluating  the  sign  of  the  determinant: 


where 
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and  X2,  l2)  and  W2  have  similar  expressions  obtained  substituting  in  the  above  determi¬ 
nants  03  with  04,  63  with  64  and  C3  with  C4.  Evaluating  the  sign  of  A  is  equivalent  to 
evaluating  the  signs  of  Wi,  W2  and  of  A'.  In  the  rest  of  this  proof  we  show  that  evaluating 
the  sign  of  A'  is  a  computation  with  degree  64.  By  using  the  same  technique,  one  can  easily 
see  that  evaluating  the  signs  of  Wi  and  W2  is  a  computation  with  degree  12. 
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We  have 


A'  =  XgiYiWl  -  Y1W2)  -  yg{XiW2  -  X2W1)  +  {X2Y1  -  X1Y2)  (1) 

By  using  the  notation  introduced  in  Section  2.1,  we  can  rewrite  Xi,  and  Yi  zs  o:  pi  -{• 
a^p2  +  as  a^pi  +  oi^p2  +  a^ps,  X2  and  Y2  as  a^pi  +  ot^ P2  +  oi^P4^  and  W2  as 

a^pi  +  a^p2  +  o:^P4f  where  pi  =  +  6?  {i  =  Considering  that  Xg  and  y,  are 

expressions  of  type  a,  and  applying  repeatedly  Rules  (1)  and  (2),  we  obtain  the  expression 


a®  +  Ol^p\P2  +  Ot^PlPZ  +  <^^PlP4  +  P2PZ  +  0^^  P2P4  +  Ot^P3p4- 


By  means  of  the  rewriting  rules  of  Section  2.1  we  have: 


Q®  +  a^Plp2  +  Ck^PlPZ  +  01^  pip  4  +  <^^P2P3  +  0‘^P2P4  +  PZP  4 

(q®  +  a^p2P3  +  C1^P2P4  +  ~  {oi^P\P2  +  PlPZ  +  P1P4Y 

^16  ^  a^^P2Pz  +  a^^p2P4  +  Ol^^P3p4 
(o'®  +  a^‘^P2P3Y  -  («’'^P2P4  +  a^^P3P4)^ 

^32  a^°p2P3 


,(7) 

^(1.2.3, 4) 
).(1,2,3,4) 

,(7) 


□ 


An  algorithm  for  proximity  queries  on  a  set  S  of  point  and  segment  sites  in  the  plane  is 
said  to  be  conventionalii  it  computes  the  explicit  extended  Voronoi  diagram  V'(5)  of  S  and 
then  performs  point  location  queries  on  V'{S)  with  a  native  method.  Note  that  the  class  of 
conventional  proximity  query  algorithms  includes  all  the  efficient  algorithms  presented  in 
the  literature.  A  conventional  proximity  query  algorithm  needs  to  perform  test  primitive 
left-right(y,  c).  Thus,  by  Lemma  24  we  conclude: 

Theorem  8  Conventional  algorithms  for  the  nearest  neighbor  query  problem  on  a  set  of 
point  and  segment  sites  in  the  plane  have  degree  at  least  64. 

Our  analysis  shows  that  performing  point  location  in  an  explicit  Voronoi  diagram  of 
points  and  segments  is  not  practically  feasible  due  to  the  high  degree. 


5.3  Implicit  Voronoi  Diagrams 

The  definition  of  the  implicit  Voronoi  diagram  V*(5)  of  a  set  of  5  of  point  and  segment  sites 
is  a  straightforward  extension  of  the  definition  for  Voronoi  diagrams  of  point  sites  given  in 
Section  3.3.  Namely  V*(5)  stores  the  topological  structure  of  the  extended  Voronoi  diagram 
V'{S)  of  5  (e.g.,  the  data  structure  of  [20])  and  the  following  geometric  information  for  each 
vertex  and  edge: 
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•  For  each  vertex  v  of  V'(5),  V*(5)  stores  the  semi-integer  (6-t-  l)-bit  approximations 
x'‘(v)  and  y*(u)  of  the  x-  and  y-coordinates  of  v. 

•  For  each  non-horizontal  edge  e  of  V'(S),  y*(5)  stores  the  pair  of  sites  £(e)  and  r(e) 
such  that  e  is  a  portion  of  the  bisector  of  £(e)  and  r(e),  and  £(e)  is  to  the  left  of  r(e). 

In  the  implicit  Voronoi  diagram  V*(S)  of  S,  test  left-right(g,  e)  is  implemented  by  com¬ 
paring  the  distances  of  query  point  g  from  sites  £(e)  and  r(e)  with  one  of  the  following 
tests,  dpending  on  the  type  (point  or  line)  of  sites  £(e)  and  r(e):  point-to-lines  distance  test, 
point-to-point-line  distance  test,  or  point-to-points  distance  test.  Thus,  by  Lemmas  2—4,  we 
have. 

Lemma  25  For  any  native  method  on  a  class  of  maps  that  includes  extended  Voronoi 
diagrams  of  point  and  segment  sites  in  the  plane,  a  point  location  query  in  an  implicit 
Voronoi  diagram  has  degree  6  and  has  the  same  asymptotic  time  complexity  as  a  point 
location  query  in  an  explicit  Voronoi  diagram. 

The  implicit  Voronoi  diagram  can  be  constructed  in  O(ralogn)  expected  running  time 
by  using  the  randomized  incremental  algorithm  of  [10].  The  most  demanding  operation  is 
incircle  test  for  three  segments,  which  has  degree  40  by  Lemma  7  (see  also  [8]).  It  is  not  hard 
to  show  that  both  the  y-ordering  of  the  vertices  of  V (S)  and  the  semi-integer  approximation 
of  the  vertex  coordinates  can  be  performed  without  affecting  the  computational  cost  and 
the  degree  of  the  computation  of  V(5). 

Lemma  26  The  implicit  Voronoi  diagram  of  a  set  of  n  point  and  segment  sites  in  the  plane 
can  be  computed  in  0(n  log  n)  expected  time,  0{n)  space,  and  degree  40. 

Theorem  9  Let  S  he  a  set  ofn  point  and  segment  sites  in  the  plane.  There  exists  an  0{n)- 
space  data  structure  for  S  that  can  be  computed  in  O(nlogn)  expected  time  with  degree  40 
and  supports  nearest  neighbor  queries  in  O(logn)  time  with  degree  6. 


6  Simplified  Implicit  Voronoi  Diagrams 

In  this  section,  we  describe  a  modification  of  implicit  Voronoi  diagrams  of  point  sites  that 
allows  us  to  reduce  the  degree  of  the  preprocessing  task  from  5  to  4  when  the  sites  are  in 
the  plane  (see  Theorems  3-5),  and  from  7  to  5  when  the  sites  are  in  three-dimensional  space 
(see  Theorem  7).  This  modification  also  has  a  positive  impact  on  the  space  requirement  of 
the  data  structure  and  on  the  running  time  of  point  location  queries. 

Let  V(5)  be  the  Voronoi  diagram  of  a  set  S  of  point  sites  in  the  plane.  We  recall  our 
standard  assumption  that  all  input  parameters  —  such  as  coordinates  of  sites  and  query 
points  —  are  represented  as  6-bit  integers. 

An  island  of  V{S)  is  a  connected  component  of  the  map  obtained  from  V{S)  by  removing 
all  the  vertices  with  integer  y-coordinate  and  all  the  edges  containing  a  point  with  integer 
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j/-coordinate.  Note  that  for  any  two  vertices  vi  and  V2  of  an  island,  y*{v-i)  —  y*{v2)  —  m  +  2 
for  some  integer  m,  where  y*{v)  is  the  semi-integer  approximation  defined  in  Section  3.3. 

The  simplified  implicit  Voronoi  diagram  V°(S)  of  S  is  a  representation  of  the  Voronoi 
diagram  V  (S)  of  5  that  consists  of  a  topological  component  and  of  a  geometric  component. 
The  topological  component  of  V°{S)  is  the  planar  embedding  obtained  from  V{S)  by  con¬ 
tracting  each  island  of  V (S)  into  an  alias  vertex.  The  geometric  component  of  V’°(5)  stores 
the  following  geometric  information  for  each  vertex  and  edge  of  the  embedding: 

•  For  each  vertex  v  that  is  also  a  vertex  of  V(S),  V°(5)  stores  the  (6-l-l)-bit  semi-integers 
approximations  x“‘{v)  and  y*{v). 

•  For  each  alias  vertex  a,  which  is  associated  with  an  island  of  V(5),  y°(5)  stores 
semi-integer  y*(a)  such  that  y*(a)  =  y*{v)  for  each  vertex  v  of  the  island. 

•  For  each  non-horizontal  edge  e  that  is  also  an  edge  of  y(S),  y°(S)  stores  the  pair  of 
sites  i(e)  and  r(e)  such  that  e  is  a  portion  of  the  perpendicular  bisector  of  £{€)  and 
r(e),  and  £(e)  is  to  the  left  of  r(e). 

The  space  requirement  of  the  simplified  implicit  Voronoi  diagram  is  less  than  or  equal 
to  that  of  the  implicit  Voronoi  diagram,  since  each  island  is  represented  by  a  single  alias 
vertex  storing  only  its  semi-integer  y-approximation.  We  can  show  examples  where  the 
simplified  implicit  Voronoi  diagram  of  n  point  sites  has  0{n)  fewer  vertices  and  edges  than 
the  corresponding  implicit  Voronoi  diagram. 

The  following  lemmas  extend  Lemmas  12-13  and  can  be  proved  with  similar  techniques. 

Lemma  27  For  any  native  method  on  a  class  of  maps  that  includes  monotone  maps,  a  point 
location  query  in  a  simplified  implicit  Voronoi  diagram  has  optimal  degree  2  and  executes 
a  number  of  operations  less  than  or  equal  to  a  point  location  query  in  the  corresponding 
explicit  Voronoi  diagram. 

Lemma  28  The  simplified  implicit  Voronoi  diagram  of  n  point  sites  in  the  plane  can  be 
computed  in  0{n\ogn)  time,  0{n)  space,  and  with  degree  4- 

The  main  advantage  of  the  simplified  implicit  Voronoi  diagram  with  respect  to  the  degree 
cost  measure  is  that  the  additional  test  primitive  needed  in  the  preprocessing  that  consists 
of  comparing  the  y-coordinates  of  two  Voronoi  vertices  (see  the  proof  of  Theorem  3)  is  now 
reduced  to  the  comparison  of  two  (6  +  l)-bit  semi-integers,  and  thus  has  degree  1.  Hence, 
the  preprocessing  for  point  location  using  a  native  method  for  monotone  maps  has  degree  1. 

By  the  above  discussion  and  Lemmas  27-28,  we  obtain  the  following  theorem  that 
improves  upon  Theorem  3. 

Theorem  10  Let  S  be  a  set  of  n  point  sites  in  the  plane.  There  exists  an  0{n)-space 
data  structure  for  S,  based  on  the  simplified  implicit  Voronoi  diagram  V°(5),  that  can  be 
computed  inO{n  log  n)  time  with  degree  4  and  supports  nearest  neighbor  queries  in  O  (log  n) 
time  with  optimal  degree  2. 
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Using  a  similar  approach,  we  can  define  simplified  implicit  order-A:  Voronoi  diagrams  for 
point  sites  in  the  plane  and  simplified  implicit  Voronoi  diagrams  for  point  «tes  in  thre^ 
dimensional  space.  This  reduces  the  degree  of  the  preprocessing  from  5  to  4  in  Theorems  4  5, 
and  from  7  to  5  in  Theorem  7. 


7  Further  Research  Directions 


Within  the  proposed  approach,  this  paper  only  addresses  the  issue  of  the  degree  of  test 
computations,  and  illustrates  its  impact  on  algorithmic  design  in  relation  to  a  central  prob¬ 
lem  in  computational  geometry.  However,  several  important  related  problems  need  further 

investigation,  and  will  be  reported  on  in  the  near  future.  ^ 

First  and  foremost,  since  the  degree  of  an  algorithm  expresses  worst-case  computational 
requirement  occurring  in  degenerate  or  near-degenerate  instances,  special  attention  naust 
be  devoted  to  the  development  of  a  methodology  that  reliably  computes  the  sign  ol  an 
expression  with  the  least  expenditure  of  computational  resources.  This  involves  the  use 
of  “arithmetic  filters,”  possibly  families  of  filters,  of  progressively  increasing  power  that, 
depending  upon  the  values  of  primitive  variables,  carefully  adjust  the  computational  effort 

(see,  e.g.,  [4,  10,  25,  34]). 

Next,  one  should  carefully  analyze  the  precision  adopted  in  executing  constructions,  so 
that  the  outputs  are  within  the  precision  bounds  required  by  the  application.  In  addition, 
each  construction  algorithm  should  be  accompanied  by  a  verification  algorithm,  which  not 
only  checks  for  topological  compliance  of  the  output  with  the  generic  member  of  its  class 
(as,  e.g..  a  Voronoi  diagram  must  have  the  topology  of  a  convex  map)  as  illustrated  in  [45], 
but  more  specifically  verifies  its  topological  agreement  with  the  structure  emerging  from 

the  tests  executed  by  the  algorithm  [38].  _ 

Beyond  these  general  methodological  issues,  the  investigation  reported  in  these  pages 
leaves  some  interesting  open  problems,  such  as  answering  nearest  neighbor  queries  in  sub¬ 
quadratic  time  and  optimal  degree  for  a  set  of  points  in  three-dimensional  space,  or  improv¬ 
ing  the  efficiency  of  the  preprocessing  stage  in  computing  the  implicit  Voronoi  diagram  of 


a  set  of  sites.  .  r  ... 

We  mention,  in  this  respect,  how  the  degree  can  impact  the  design  of  geometric  primi¬ 
tives  adopted  in  existing  algorithms  for  Voronoi  diagrams  of  point  and  segment  sites.  Let 
(01,02,  as;  04),  with  oi  either  a  point  or  a  segment,  denote  the  incircle  test,  where  04  is 
tested  for  intersection  with  circle(ai,a2,a3).  Specifically,  consider  (pi,P2, hiPs),  which  can 
be  answered  with  degree  12  [8].  We  show  that  it  can  be  more  efficiently  executed  as  fol¬ 
lows.  Perform  first  the  test  (pi,P2,P3; h)-  Let  c  and  c'  be  the  centers  of  circle(pi,p2, ^i) 
and  circle (pi,p2, Pa) >  respectively.  Two  cases  are  possible;  either  circle (pi,p2,P3)  intersects 
li  or  it  does  not.  In  the  first  case,  p3  is  inside  circle(pi,p2,Ii)  if  and  only  if  c'  and  ps  lie  on 
opposite  sides  of  fine  through  pi  and  p2  (see  Figure  3  (a)  and  (b)).  In  the  second  case 
the  answer  to  test  (pi,P2,  IilPa)  depends  on  which  side  of  piP2  point  pz  lies  (see  Figure  3 
(c)  and  (d)).  Thus  test  (pi,P2,  IiiPa)  is  reduced  to  test  {puP2,P3\h)  that  can  be  executed 
with  degree  8  (see  [8]),  and  at  most  two  other  left-right  tests  of  lower  degree. 
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(c)  (d) 

Figure  3:  Different  cases  for  test  ipi,P2th'iP3)- 


Finally,  an  important  issue  for  future  research  deals  with  the  experimental  comparison 
between  point  location  algorithms  in  implicit  Voronoi  diagrams  and  traditional  point  loca¬ 
tion  algorithms  in  explicit  Voronoi  diagrams.  We  are  currently  implementing  GeomLib,  an 
object-oriented  library  for  robust  geometric  computing  that  will  be  accessible  through  the 
World  Wide  Web  by  using  the  architectural  framework  of  Mocha  [6,  5]. 
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